Analytic Solution of Strongly Coupling Schrodinger Equations 
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A recently developed expansion method for analytically solving the ground states of strongly 
coupling Schrodinger equations by Friedberg, Lee and Zhao is extended to excited states and applied 
to the pedagogically important problems of power-law central forces. With the extended method, the 
Hydrogen atom problem is resolved and the low-lying states of Yukawa potential are approximately 
obtained. 
PACS: 03.65. Ge 



oo 



o 
S3 

Oh 



=5 



I. INTRODUCTION 
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(N ! The main task in application of quantum mechanics is to solve Schrodinger equations for various potentials. Un- 
fortunately, only few of them can be solved exactly. In textbooks various approximation methods, especially the 
perturbation theory, are discussed with great length. While the perturbation theory, as the most important approx- 
imation method, has exhibited its power in a lot of weakly coupling cases, it is not applicable to strongly coupling 
problems which, however, exist widely and play important roles in the real world. How to analytically solve the 
Schroedinger equations in nonperturbative case remains an open question. 

Recently, a new method to discuss strongly coupling potentials was presented by Friedberg, Lee and Zhao [1] . The 
goal is to solve the ground state of the time independent Schroedinger equation with strongly coupling potential 

> : 

rn '■ + v)i> = Ei, . (i) 

Considering a boundary potential V(f) > with the minimum V(r) — at the origin point r — 0, one introduces a 
| dimensionless scalar parameter g which reflects the intensity of the potential by writing 

V(f) = g 2 v{r) . (2) 



For strong potentials we are interested in, g is a large number, the ground state ip and the corresponding ground 
energy E can be expanded in the inverse number of g, 



S(r) = gS (f) + <?%(f) + g- 1 S 2 (f) + ■■■ 

E = gEa+g a E 1 +g- 1 E 2 + --- . (3) 



Substituting the expansion into (1) and then comparing the coefficients of g "on the both sides, one obtains a series 



of first order differential equations (by contrast, the Schrodinger equation itself is of second order): 

(yS ) 2 = 2mv, 
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X/Sq ■ ySi = - v Sb - m£; o, 
VSo • yS 2 = t;(V 2 Si - ySi • VSi) - mEi, 



(4) 



where m is the particle mass. The first equation of the above hierarchy contains only So as unknown term and can 
be solved by direct integration with boundary condition. The second equation contains only Si as unknown term 
(note that So has been solved by the first equation) and can be solved in the same way as So, and E can be fixed 
by considering the convergence of Si at the origin. With the same procedure, all Sj and Ei can be determined 
step by step. To illustrate the procedure clearly, one considers the familiar one dimension harmonic oscillator with 
V(x) — g 2 m 3 x 2 . In this case the hierarchy (4) is reduced to 
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Q 2 = 2mV, 
dS dSi 1 d 2 S 



dx 2 dx 2 

dS d5 2 l,d 2 S*i 





dx dx 2 dx 2 



(5) 



Using the boundary condition tp(±oo) = or 5o(±oo) = — oo and the convenient renormalization ^>(0) = 1 or 
5o(0) = 0, the integration of the first equation gives So(x) = — -^m 2 x 2 . With the solved So, the second equation is 

reduced to S[ = (| — -j^)/ x - To avoid divergence of S[ at x = 0, the only way is to choose E — -j=- With the same 
procedure, all Si and Ei with i > can be proven to be zero, and the final result gives precisely the same solution as 
that in any textbook. For details, see [1]- [3]. 

The readers who are familiar with the well-known WKB approximation may ask the difference between the WKB 
and the hierarchy (4). While the hierarchy (4) looks similar to the WKB hierarchy [4], they are very different from 
each other. The most explicit difference is that in WKB the energy E appears always in the fist equation of the 
hierarchy together with the potential V, but here E and v may not appear in the first equation simultaneously, see 
(4). In WKB, the energy is determined by the continuity condition of wave functions at the so-called "classical turning 
point" , or equivalently by using the Bohr-Sommcrfcld like quantization condition. In the approach here, the energy 
eigenvalue is directly determined by the boundary condition. Physically, in WKB the wave function is expanded in 
U, it is useful only for semiclassical systems, while the wave function and energy eigenvalue here are expanded in 
the inverse number of the strong coupling constant, the hierarchy (4) is like the expansion in the inverse number of 
large N often used in nonperturpative treatment of quantum field systems, see, for instance, [5]. Also the WKB wave 
function e^ a ^ is an oscillating form which corresponds to real particles, so that more applicable to semiclassical case, 
but the wave function in (4) is a damping form and more suitable for dealing with bound states. 

Although the new method is powerful in solving strongly coupling Schrodingcr equations, the step of writing the 
wave function ijj(r) in the form ?/>() = e s ^ in (3) is a strong constraint in the application of the method to more kinds 
of potentials. For example, the linear potential V — gr cannot be solved through expansion (3) due to the divergence 
encountered in dealing with (4) . From the analysis in Section 2 we will see that the expansion with the form of (3) is 
only suitable for harmonic-oscillator-like potentials which can be expanded quadratically around the minimum. On 
the other hand, the expansion (3) is valid mainly for ground states, an excited state can not be expressed as a pure 
exponential function. 

To make the expansion method more applicable, one should construct carefully the expansion form of tp(r) through 
physical analysis of different potentials and of different states. In this paper, we first investigate the coupling constant 
dependence of bound states for power-law central forces in order to construct proper expansions of ip and E in 
inverse coupling constant in Sections 2 and 3, and then apply the extended expansion method to Coulomb potential 
and resolve the excited states in Section 4. Finally we discuss Yukawa potential and obtain the low-lying states 
approximately in Section 5. The Conclusions are given in Section 6. We have chosen K = c = 1 throughout the paper. 



As mentioned above, the dependence of bound states solutions on coupling constant is crucial since the method 
is based on the coupling constant expansion. Qualitatively speaking, when coupling constant g in potentials like (1) 
increases, the bound state energy drops down and the space extension of the corresponding wave function becomes 
narrower, that is, the system is more bounded. In the following we concentrate our discussion on power-law central 
forces V(r) oc g k r n , the reason is twofold: first, with those potentials the coupling constant dependence can be 
factorized, as we will see below; second, any boundary potential closely around its minimum point can be approximated 
well by power-law curve with certain value of n. 



II. COUPLING CONSTANT DEPENDENCE OF BOUND STATES SOLUTIONS 



A. Coupling constant dependence of energy 



Let's consider a power-law central potential 
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- p 2 

2m 

V(r) = ±g k m(mr) n . (6) 

In order to keep boundary condition we should choose the sign + for n > and — for n < 0. Here the coupling 
constant g is guaranteed to be dimensionless. For instance, the harmonic oscillator V — \mbj 2 r 2 can be written as 
V = g k m 3 r 2 with g k = \{^) 2 ■ 

Considering g as a parameter in the Hamiltonian, we apply Hellmann-Fcynman Theorem [4] for any arbitrary 
energy level E and the corresponding state \ip) 

to (6) and obtain 

dE k 



From Virial Theorem [4], 



<hi g (wm- (8) 



we get 

dE 2k 1 



E, (10) 



dg n + 2 g 
and its solution 

E(g)=g#*e, (11) 

where e depends on r only, the coupling constant dependence is factorized. A useful consequence of this factorization 
is that for a given potential V(r) one can determine firstly e in some ideal case, for instance, the limit g 3> 1 or g <C 1, 

2k 

and then obtain the energy E of the real system by multiplying e by the g-factor g™+ 2 . 

When a potential is not in the power-law form, the factorization (11) fails. However, in strongly coupling cases the 
low-lying states are restricted in a small region around the origin r = 0, and the potential can be expressed in the 
power-law form 

lim V(r) oc ±g k m(mr) n (12) 

in the neighborhood of the origin, the factorization is still approximately valid, and the leading term E in the 
expansion (3) is given by (11). This conclusion will be used to solve the Yukawa potential 

V(r) = -g 2 e- ar /r (13) 

in Section 5. 

B. Scale transformation 

For central forces, the bound state wave functions can be written as [4] 

il>(f>) = R(r)Y L M(6,4>), (14) 

where Ylm(9,(/)) is the spherical harmonic function carrying quantum numbers L and M, and R(r) which is related 
to the coupling constant satisfies the radial equation 

d 2 R(r,g) 2dR(r,g) „ , „ T , N , L(L + 1) „, N „ 

dr 2 ^ + rdT^ + 2m{E V)R[r ' 9) — - R ( r > 9) = 0. (15) 
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Considering a scale transformation 



r — > r s = ar, 

9 -»• 9s = bg, (16) 

with a and b being two positive and real constants, and using the factorization (11) the radial equation (15) is 
transformed into 

J_^r^) + J_2 dRir^ + ^ _ , _ J_X(Z^1) = Q _ 

If the constants a and 6 are restricted by 

b^h =b k a n =a- 2 , (18) 

R(r S7 g s ) and i?(r, g) satisfy the same radial equation and the same boundary conditions at the origin and at infinity ( 
the transformation does not change boundary conditions ), that is 

R(r s ,g s ) = cR(r,g), (19) 

where c is a g and r independent constant and can be removed by normalization. 

The invariance of the equation of motion under the scale transformation (16) and (18) originates from the symmetry 
of the Hamiltonian (6), 

H^H S = \A. (20) 
a 1 

Under the transformation all eigenvalues of the Hamiltonian are rescaled by a constant factor 1 /a 2 and all eigenfunc- 
tions are unchanged. 

The scale transformation property of the radial function implies that we can replace the original two variables r 
and g by one dimensionless scale variable 

s = g k (mr) n+2 (21) 

which is invariant under the scale transformation (16) and (18). With properly selected normalization scheme, the 
radial functions R(si) and R(s2) with two different coupling constants g\ and g<i are related by 

R{g k {mr) n + 2 ) = R (g k (|) (mr)»+ 2 





(22) 

Therefore, the wave function at point r with coupling constant g\ equals the wave function at point {gi/g2) k ^ n+2 r 
with 52 • The wave function is contracted in a smaller region when the coupling becomes stronger while extended to 
a larger region when the coupling becomes weaker. With the relation (22), the behavior of a weakly coupling system 
at the place far from the origin is equivalent to that of a strongly coupling system near the origin. 
The conclusion of this section can be summarized as 

V = ±g k m(mr) n , 

2k 

E = gn+2 6 , 

R(r,g) = R(g k (mr) n+2 ). (23) 

For the potentials with power-law behavior (12) around the origin, the leading terms of the energy and wave function 
are also given by (23). This conclusion tells us how to write down the expansion in g for a power-law or a power- 
law-like central potential. For example, the expansion (3) is only suitable for harmonic oscillator and the similar 
potentials. 
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III. EXTENDED EXPANSION METHOD 



The key point of the method by Friedberg, Lee and Zhao is to expand the energy and wave function of a strongly 
coupling system in the inverse number of the coupling constant g. As mentioned above, the original expansion (3) is 1) 
only suitable for the potentials with harmonic oscillator behavior around the minimum, and 2) only valid for ground 
states. To solve (1) by the expansion method, we should proceed three steps: First, write a proper g-power expansion 
for energy and wave function according to the scaling property discussed in previous section; Second, substitute the 
expansion back into (1) and obtain a hierarchy of first-order differential equations; Last, integrate the hierarchy one 
by one with the help of the boundary and convergence conditions, as described in the example of one dimension 
harmonic oscillator in Section 1. 

As a pedagogical goal, let's extend the method to solve three dimensional Coulomb and harmonic oscillator problems 
in quantum mechanics. 

In [1] and [2], the expansion and solution of the ground state are given by 

v-S- 

r 

E = g 4 E a + g 2 E 1 + g°E 2 + ■■■ 

1 4 

= q m 

2 y 

S(r) = g 2 S (r)+g°S 1 (r)+g- 2 S 2 (r) +■■■ 

= g 2 mr (24) 



for Coulomb potential, and 



V= l q 2 m 3 r 2 

2 y 



E = gE Q +g°E 1 +g- 1 E 2 + --- 
3 

= 2 9m 

S(r) = gS (r) + <;%(r) + g' 1 S 2 (r) + ■ ■ ■ 

= \gm 2 r 2 (25) 

for harmonic oscillator. It is clear that the expansions for the two potentials are in accord with the scaling law shown 
in (23). 

Since any excited state is not so bounded like the ground state, we assume that the excited states differ from the 
ground state e^ s ^ by a smooth function P(r), 

R(r) = P{r)e- S ^ . (26) 

When P(r) = 1 we go back to the ground state. From the spirit of the expansion method we express P(r) and S(r) 
as finite polynomials of l/g. With the help of the scaling law (23), they can be further rewritten as polynomials of a 
dimcnsionless scaling variable g 2 mr for Coulomb potential 



i 2 mrY 



or gm 2 r 2 for harmonic oscillator 



-7 

1=0 

E - g 4 e , (27) 
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E = ge. (28) 

Here a,/3, 7 are nonnegative integers, and dj,&i_j are g and r independent constants. Note that the expressions (27) 
and (28) obtained from the scaling analysis guarantee a wider spread of excited states than the ground state. The 
expansions here are also in accordance with the discussion with algebra method [6] . 

We now substitute the expansions into the corresponding radial equations and discuss the convergence in the limit 
r — > 0. It is easy to prove that when a > the most divergent term in the radial equation is a 2 ] g~ 4a b 2 c+1 r~^ +2a+2 '> 
for Coulomb potential and 4a 2 g~ 2a b 2 l:+1 r~ 2 (' 1+2a+1 * > for harmonic oscillator. Since there is only one most divergent 
term, it is impossible to cancel it by other terms. The only way is to take a = 0. In this case, the expansion of S(r) 
in (27) or (28) is reduced to two terms, and the constant term can be absorbed by the normalization of the wave 
function. Moreover, the terms with negative powers (7 > 0) in the expansion of P(r) must vanish also in order to 
keep the whole radial wave function R(r) — P(r)e~ s ^ from divergence at the origin r — 0. Finally we obtain for any 
bound state the extended expansion 

S(r) = b g 2 mr , 

N-l 

P N(r) = a i(9 2 mrY (29) 

i=0 

for Coulomb potential, and 

S(r) = b gm 2 r 2 

N-l 

P N (r) = ]T a k (gm 2 r 2 y (30) 

i=0 

for harmonic oscillator with N > 1. 



IV. SOLUTIONS FOR COULOMB POTENTIAL 



We now turn to the second step, namely substitute the convergent expansion (29) into the radial equation for 
Coulomb potential, and get a series of equations in different orders of g: 



-2AT+2 



,27V 



g 2k (0 <k<N) 



b\m + 2e = 
1 - b N = 

{k(k + 1) - L(L + 1)) a k + ^(N - fc)o fc _i = 

L{L + l)a = 



(31) 



From the first two equations we derive immediately the energy level e^v = — j^i and the coefficient bo = jj in S. 

To obtain the coefficients en in the polynomial P, we need to discuss the relation between the two quantum numbers 
7V(> 1) and L(> 0), they come from the polynomial Pn{t) and the spherical harmonic function Ylm(&, <ft)- For N = 1, 
the requirement Pi(r) = ao ^ leads to L = from the last equation of (31), the general expression (29) is then 
reduced to the ground state solution (24). 

For N > 1, we discuss three cases for the quantum number L separately: 1) L = 0, the third equation of (31) 
determines the recursion relation between the coefficients a k 



R 



NO 



/N-l \ 



a k 



2{N - k) 
' Nk(k+ 1) 



< k < N, 



(32) 



the only unknown coefficient a ^ is determined by the normalization. 2) < L < N, the last two equations of (31) 
lead to the coefficients 00 = ai = 02 = ■ ■ ■ = aL-i = 0, and then the radial solutions are written as 



G 



'N-l 



\k=L / 

2 N — k 

ak= N L(L + l)-k(k + l) ak ^ L<k<N (33) 

cil 7^ is determined by the normalization too. 3) L > N, the last two equations of (31) require do = a\ = a 2 = ■ ■ ■ = 
ax -i = 0, and there is no nonzero solution in this case. 

We summarize the solution of the bound states for Coulomb potential: 



1 

2N 2 

- 2 N ~ k 
ak ~ N L(L + I) - k(k + l)^ 1 ' 

N>0, L = 0,1,2,- • -,N- 1 (34) 

with a i, 7^ determined by the normalization of the wave function. (34) is exactly the same as what obtained by 
solving the second-order Schrodinger equation directly in normal textbooks [4], but here by using the expansion in 
coupling constant and with the help of the scaling law the complicated second-order differential equations are replaced 
by simple algebra equations. 

In a similar way one can also derive all the bound states for harmonic oscillator by substituting the convergent 
expansion (30) into the radial equation. 



V. YUKAWA POTENTIAL 



The potential (13) was first introduced into physics in 1930s by Yukawa in the study of strong interaction between 
nucleons through meson exchange [7] . It is also known as Debye-Hiickel potential in plasma physics and Thomas- Fermi 
potential in solid-state physics [8]. While the Yukawa potential is important in physics, the corresponding Schrodinger 
equation can not be solved analytically and exactly. A lot of publications have been contributed to the problem using 
different approximations (see [9] - [12] and references therein). Here we use the extended expansion method to obtain 
its ground state and low-lying excited states. Comparison with strict but numerical calculation will be made. 

We first consider the proper 1/g expansion for Yukawa potential. Since Yukawa potential is not a power-law central 
force, the g dependence of its bound states can not be simply factorized like (23). In strong coupling case with 
g » 1, however, the low-lying wave functions are mainly distributed in the vicinity of the origin. Since in the limit 
of r — > 0, Yukawa potential approaches to Coulomb potential, we can thus write the expansion of the low-lying states 
of Yukawa potential by properly modifying the expansion of Coulomb potential. The inverse of the parameter a in 
the Yukawa potential (13) represents the mean potential range. To make our calculation based on Coulomb potcntial( 
whose potential range approaches infinity) more effective, we should require the mean potential range to be much 
larger than the "Bohr radius" of Coulomb potential, that is, ^ >> or equivalently -^§^ « 1. We will show 

below how the dimensionless number -^-^ controls the expanding series. 

We now deal with the ground state. In the light of (24) for Coulomb potential, we expand the energy and wave 
function as 

E = g^Eo + g 2 E 1 + g a E 2 + ■■■ 
R = e" s W 

S(r) =g 2 S Q (r) +g a S 1 (r) +g- 2 S 2 (r) + ■■■ (35) 

Note that the expansions of E and S are finite for strict Coulomb potential but infinite for Coulomb-like potentials. 

By Substituting the expansions (35) into the corresponding stationary Schrodinger equation, we obtain an infinity 
hierarchy of first-order differential equations (5). With similar procedure in dealing (5), we can solve them one by 
one to any order we want. The energy obtained to 0{g~ A ) and the wave function to 0(g~ 2 ) read: 
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£ = "> + 5 2 a - -p- + .<r 2 ^ + O <T 4 = 3 4 ™ - o + — - 7 — + o — + ° — 
2 4m 2m/ 2 <?^m 4 g^m. 2 g^m g^m 



g 2 mr + (l- e — )- a]+ 0(.g- 2 )=.g 2 m[r+— / dr'(— — — - 1) + 0{( -^-) 2 )]. (36) 

Jo ' 9 171 Jo ar g m 



We see clearly that the dimensionless parameter a/g 2 m controls the degree of convergence of the expansion. When a 
vanishes, Yukawa potential is reduced to Coulomb potential, and therefore, the above solution is reduced to the strict 
solution (24). 

We now consider the lowest excited state. By modifying the expansion (29) with N = 2 for Coulomb potential, we 
write the following expansion 

E = g A E + g 2 E 1 + g°E 2 + ■■■ 
yj = P(r)e- s ^Y LM (d,0) 
P(r) = g 2 h(r) + g%(r) 

S(r) = g 2 So(r)+g a S 1 (r)+g- 2 S 2 (r) + ■■ , (37) 

Unlike the expansion for strict Coulomb potential, here E and S{r) have infinite terms, and the r-dependence of &o 
and b\ in P(r) is unknown and only in the limit case we have: 

lim bi(r) = lim bAr) oc mr , 

lim b (r) = lim b (r) oc 1 . (38) 

Substituting the expansions (37) into the radial equation with N = 2, we get again a series of first-order differential 
equations. The first two which are respectively proportional to g 6 and g read 

dS 



dr = ^-2mE , 
dS dSi IdbidSo 1 dS _ 

— — — = —- — _ + -(_ me ar )-mE 1 . (39) 

dr dr b\ dr dr r dr 

By using (38) the condition to keep convergence of dSo /dr and dS\ /dr at r — result in 

Ef\ = , 

8 

dSo m 
~dr~ ~ ~2 ' 

f^I^fZl + J.* 2£l . (40) 
dr r bi dr 

In a similar way we can solve the equations proportional to g 2 and g° to determine E\ and S\ 

_ 6 + L(L+l) 
Ei - a , 

6 = 1 - L , 

61 = -^(l-e— ). (41) 

Put the constituents obtained above together we finally get approximately the analytical solution of the lowest excited 
state for Yukawa potential as follows 



E = + 6 + L [ L + 1) g 2 * + Otf) = -^[1 - -£- ( I, 2 + L + 6) + 0((^-) 2 )] 

r rv f r 1 — 9p- ar ' p- a ' 

S(r)=g 2 m[- + ^ dr' C , + 

2 g z m J ar' 1 — e 



g 2 m g 2 m 



2^ 

.g z m' 



P(r) = a ((l-£)-^(l-e— ) 
N = 2, L = 0,1, (42) 
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the only constant ag is determined by normalization. Again, the dimensionless parameter a/g 2 m dominates the 
convergence. The analytic result E for ground and lowest excited states is compared with the strict but numerical 
result En and with the Coulomb result Ec in Tabel (1) for different values of the dimensionless parameter a/g 2 m. 
All the energy values have been scaled by g 4 m. Since the ground state energy is calculated to the fourth order, it 
always agrees well with the strict solution. The Coulomb energy is just the leading order of the Yukawa energy, 
its deviation from the strict result is very strong, especially for the excited states, then relatively large. When the 
parameter a/g 2 m decreases, the analytic result to the second order for the excited states becomes more and more 
close to the strict solution. 

We also compared our energy eigenvalues with the often cited numerical calculations in [12]. The deviation of our 
analytical result from the high accurate numerical result is less than 2.5% for Is state with < 1/3, 5% for 2s state 
with < 1/30, and 5% for 2p state with < 1/20. The deviation decreases when the dimensionless parameter 
a/g 2 m becomes smaller. 

The other excited states ( N = 3, 4, • • • ) can be obtained with the similar method. 
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FIG. 1. Comparison of our analytic result E of Yukawa potential to the fourth order for the ground state and to the second 
order for the lowest excited states with the strict but numerical result En and with the Coulomb result Ec for three values of 
the dimensionless parameter a/g 2 m. 

The solutions (36) and (42) have the following features: 

1) The leading orders are just the solutions of Coulomb potential, and therefore are invariant under the scale trans- 
formation (16). The high order corrections break the scale invariance explicitly. When a — > or r — > 0, the high 
order corrections approach to zero and the scale invariance is restored. 

2) Different from the accidental degeneracy for Coulomb potential, the energy level for Yukawa potential depends on 
the quantum number L. 

3) The expanding series are controlled by the combined dimensionless parameter For sufficiently small values, 
namely for very strong interaction, the expansions (36) and (42) converge very fast, one can consider the lower-order 
contributions only. 

VI. SUMMARY 

The recently developed expansion method [1] - [3] provides an alternative way to solve analytically Schroedinger 
equations with strong coupling. The key point of the method is to construct proper expansions for different potentials 
and different quantum states. Through investigating the scale transformation invariance of the bound states for 
power-law central forces, we have extended the expansion method from ground state to any excited state for power- 
law central forces and to low-lying states for powcr-law-like central forces. With the extended expansion method, we 
obtained analytically the strict solutions of all bound states for Coulomb potential and the approximate solutions of 
low-lying states for Yukawa potential. Further application of the method to other physical potentials is of interest 
and value. 
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